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Abstract 

Many electrical grid transients can be described by 
the propagation of electromechanical (EM) waves that 
couple oscillations of power flows over transmission 
lines and the inertia of synchronous generators. These 
EM waves can take several forms: large-scale stand- 
ing waves forming inter-area modes, localized oscil- 
lations of single or multi-machine modes, or traveling 
waves that spread quasi-circularly from major grid dis- 
turbances. The propagation speed and damping of these 
EM waves are potentially a powerful tool for assess- 
ing grid stability, e.g. small signal or rotor angle sta- 
bility, however, EM wave properties have been mostly 
extracted from post-event analysis of major grid distur- 
bances. Using a small set of data from the FNET sen- 
sor network, we show how the spatially resolved Green 's 
function for EM wave propagation can be extracted from 
ambient frequency noise without the need for a major 
disturbance. If applied to an entire interconnection, an 
EM-wave Green's function map will enable a model- 
independent method of predicting the propagation of 
grid disturbances and assessing stability. 



1 Introduction 

The frequency of the electrical grid is determined to 
a great degree by the rotational speed of large syn- 
chronous generators at centralized power stations!!]. 
In equilibrium, the mechanical power supplied to the 
shafts of these generators is in balance with the elec- 
trical power withdrawn by the electrical loads and sys- 
tem dissipation. Fluctuations in this power balance both 
large and small do occur, and the kinetic energy stored 
in the rotational inertia of these large generators and as- 
sociated turbines are the initial buffer against major fre- 
quency changes, e.g. if the electrical load were to sud- 



denly increase without a corresponding increase of input 
mechanical power, the power to supply the load is ex- 
tracted from the rotational kinetic energy of the genera- 
tors, causing them to slow and the electrical frequency 
to decrease. Such power imbalances are typically not 
sustained because feedback control systems, i.e. gover- 
nors and automatic generation control, modulate the me- 
chanical input power and restore balance within a few 
seconds to tens of seconds. In this manuscript, we are 
interested in time scales that overlap with both the "in- 
ertial" dynamics of the electrical system as well as the 
primary governor response!]]. 

The electrical transmission grid is a distributed sys- 
tem, i.e. the electrical generation and load are point 
sources spread over a wide geographical area (~ 
2000 km) interconnected by a grid of high voltage (> 
115 kV) transmission lines. On longer time scales, the 
frequency across this distributed system is uniform, and 
any power imbalance is initially a local phenomenon. 
When the electrical load suddenly increases at a partic- 
ular transmission bus, the phases (and power flows) be- 
tween that bus and its nearest neighbor buses suddenly 
increase. The increased power flow causes the rotational 
speed of nearby generators to slow as they deliver their 
stored kinetic energy to the increased load. This slowing 
increases the electrical phase angle between the nearby 
generators and their nearest-neighbor generators. The 
transfer of power from the nearest neighbors supports 
the original generators, but their rotational speed also 
slows (to a lesser degree as there are more of them) 
building up phase difference and power flow from more 
distant generators. The effect is an approximately circu- 
lar, outward propagating electromechanical (EM) trav- 
eling wave whose speed ( 1000 km/sec) is determined 
to first order by the rotational inertia of the generators 
and the susceptance of the transmission lines|2|. These 
properties can show significant spatial variation over an 
electrical grid, particularly in situations where genera- 
tors and loads are clustered yet separated from one an- 
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other. There can also be temporal variation on many 
time scales as patterns of generation and load change in 
time. In addition, the wave propagation may appear to 
be one or two-dimensional depending on the topology 
of the electrical grid O - 

Many transient grid phenomena can be understood 
from the standpoint of EM waves, e.g, inter-area oscil- 
lations (i.e. global modes) are a manifestation of an 
interconnection-scale standing EM wave. From typi- 
cal EM wave speeds and interconnection size, we esti- 
mate the frequency of the fundamental, one-wavelength, 
north-south mode in WECC to be approximately fys ~ 
[1000 km/sec]/[3000 km] ~ 0.3 Hz, consistent with ob- 
served values[ 1 1- The WECC modes with shorter wave- 
lengths have higher frequency] 1|, which is also con- 
sistent with the interpretation in terms of standing EM 
waves. The physics of local plant modes is the same 
as inter-area oscillations!]]], however, local plant modes 
typically only extend over a single or a few genera- 
tors and a small part of the transmission system, inval- 
idating a formal continuum description in terms of EM 
waves. In spite of this, the technique we develop in this 
manuscript should also be applicable to characterizing 
local plant modes. 

Inter-area and local plant modes are small signal os- 
cillations that, under the right conditions, may become 
unstable and spontaneously grow in amplitude. Even 
if stable, these modes can be excited to significant am- 
plitude by a major system disturbance. In either case, 
the resulting oscillations can lead to protective relaying 
actions that may cause loss of load, system separation, 
or large-scale blackoutsHHU. A good understanding of 
EM-wave propagation properties will create a better un- 
derstanding of the oscillatory modes and transients of an 
electrical interconnection. Even better, a real-time wide- 
area measurement system (WAMS) that extracts the EM- 
wave properties would provide a model-independent 
method for estimating oscillatory mode shapes, frequen- 
cies, and damping and a method for predicting how a 
major disturbance will impact the rest of the intercon- 
nection. The intent of this manuscript is to demonstrate 
the feasibility of extracting the Green's functionQ for 
EM wave propagation from ambient frequency noise, a 
crucial first step towards realizing such a real-time mon- 
itoring tool. 

2 Model-based stability 

Rigorous one and two dimensional continuum models of 
transmission and generation J21 [3] have been developed 
that describe EM waves, and they predict propagation 
characteristics similar to that described above. These 



models are useful for developing qualitative understand- 
ing of many transient grid phenomena, however, the use 
of these reduced models in grid operations and moni- 
toring is in doubt because the quantitative accuracy of 
these models is questionable due to the complexity of 
real electrical loads and control systems and the pres- 
ence of additional equipment not included in the mod- 
els, e.g. DC lines/ties, phase shifting transformers, and 
other FACTS devices. 

Specialized techniques have been developed for com- 
puting the mode frequencies and damping for detailed 
models of large interconnections with many generators 
and transmission lines|[TJ. For a given system condi- 
tion, these techniques can accurately characterize the in- 
terconnection's oscillatory modes. Also, detailed time- 
harmonic transient models of an interconnection can 
capture the details of EM waves, i.e. a frequency distur- 
bance propagating outward from the initiating event|6l 
even including the reflections of the waves off the grid 
boundaries. We note that the accuracy of both types 
of off-line studies is dependent on the time-consuming 
task of creating models with accurate parameters for 
many tens of thousands pieces of grid equipment, how- 
ever, the uncertainty in model parameters due to con- 
tinually changing grid operating conditions and equip- 
ment availability calls into question the viability of real- 
time assessment of interconnection stability using these 
model-based techniques. Here, we explore the feasibil- 
ity of estimating EM-wave Green's functions from am- 
bient frequency noise in WAMS measurements. These 
Green's functions will form the basis of a future model- 
independent method for assessing small signal and tran- 
sient stability. 



3 Existing WAMS measurements 

The advent WAMS[[7][E) allows for detection and mon- 
itoring of EM waves in much more detail than in the 
past. WAMS-based study of EM waves has focused 
on two major areas: the propagation and spatial lo- 
calization of major frequency disturbances (9] [10] [7] |6l 
and the damping and mode structure of major inter-area 
oscillations ifTTI fl2l [T3l . In the first example, a sudden 
loss several hundred MWs of generation due to a sys- 
tem disruption results in a frequency decline of 0. 1 Hz, 
which is easily detectable by WAMS lfTOl . Measuring 
the differences in arrival times at different points in the 
WAMS and assuming or estimating EM wave propaga- 
tion speed allows for approximate location of the ini- 
tial disturbance[9|. Alternatively, the WAMS transient 
data could be used to estimate EM wave speeds between 
the known disturbance location and the detection point, 
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however, major frequency disturbances only occur once 
every few days and the conditions of the electrical grid 
are continually changing, which makes the applicabil- 
ity of WAMS -extracted EM wave properties uncertain 
for grid conditions as little as an hour later. In the fu- 
ture, this time scale will likely shorten to a few min- 
utes as time-intermittent renewable generation, such as 
wind and solar photovoltaic, will create fast and stochas- 
tic changes in grid conditions. It is exactly under these 
highly variable conditions that real-time monitoring of 
grid conditions via EM waves will be the most valuable. 

In the second example, correlations in the frequency 
spectra of noise-driven (or probe-signal driven) fre- 
quency oscillations are used to identify the amplitude 
and time-phase of standing EM-wave oscillations at sys- 
tem buses participating in a small number of inter-area 
oscillations 1 13 1. We note that in |[T3l . the frequency of 
the single mode that was studied was previously iden- 
tified by an eigenvalue analysis. The reliance on corre- 
lations of spectra between buses appears to restrict this 
analysis to a small numbers of modes that show a signifi- 
cant level correlation across the interconnection. Modes 
that are spatially localized will not show such long-range 
correlations and would be difficult to identify. By only 
studying a single or a few modes, this analysis ignores 
a wealth of information in the remainder of the spectra 
that describes how transients propagate throughout the 
system. 

4 Green's function estimation 

As an EM-wave propagates from a bus at location Xi 
to a bus at location Xz, it effectively encodes all of the 
grid properties between X] and X2 making EM waves 
a potentially powerful grid diagnostic. For a Dirac 
delta-function frequency disturbance at X2 at t = 0, 
i.e. A/(x,f) = A/q8(x — Xz)8(t), the time response 
of the frequency at xi is called the Green's function 
Ga/(xi,X2,?). In this manuscript, we demonstrate the 
feasibility of estimating Ga/(xi,X2,£) for EM wave 
propagation in near real time from ambient frequency 
noise. We note that the frequency disturbance used to 
define Ga/(xi,X2,?) also corresponds to a step change 
in phase at X2, i.e. A4>(x,f) = A4>o§(x — X2)®(f), which 
corresponds to an abrupt change in the local generation- 
load balance as described above and investigated in 

Our Green's function estimation approach is adapted 
from and closely follows recent theoretical work by 
Snieder[14 ] on passive imaging via correlation of earth- 
quake coda waves in locally isotropic media and re- 
lated experimental work by Derode et q/ fl5l . In this 



manuscript, we review the relevant parts of Snieder's 
work and point out the essential differences with the 
present problem. Additional references on the founda- 
tions of this method can be found in the citations of lfl4ll 
and |Q3]. 

Snieder considers estimation of the Green's function 
for wave propagation between two observation points xi 
and X2 in a two-dimensional wave propagation medium 
that is embedded with many randomly placed scatterers 
(indexed by s). Snieder assumes that a wave is launched 
into the medium by some external source or event, and 
the scatterers act as secondary sources of singly and 
multiply scattered waves such that each scatterer re- 
emits a wave given by S s (t). The S s (t) are uncorrected 
because the scattering sites are randomly placed. These 
waves propagate from their scattering (source) locations 
to xi and X2 (with their individual time delays), and the 
linear superposition of these waves form the aggregate 
signals received at xi and X2, i.e. p Xl (t) and /? X2 ( f )- 
Next, we form the cross correlation of the two received 
signals over a time window T 

C(xi,x 2) t)= [ p X2 (x+t)p Xl (x)dx. (1) 
Jo 

In principle, C(xi,X2,f) involves a double sum 
over scatters s for p Xl (t) and s' for p X2 {t), however, 
Snieder| 14 1 shows that if the individual scattered waves 
(i.e. sources) S s (t) do not have a time average compo- 
nent, the cross terms with s 7^ s' average to zero after a 
sufficiently long T, leaving only the diagonal terms with 
s = s' . This requirement on the S s (t) will be crucial in 
determining the final form of Ga/(xi,X2,?). 

With only the diagonal terms (s = s ! ) remaining, 
C(xi,X2,?) is simply a sum of autocorrelations of the 
individual source signals S s (t) where the time delay in 
the autocorrelation is a sum of the original delay t in the 
cross correlation from Eq.[T]and the difference in the ar- 
rival times at xi and X2 from scatterer s. The fact that 
the S s (t) do not have a DC component means that they 
are oscillatory in nature, therefore, the individual auto- 
correlations will be oscillatory are well. The phase of 
these oscillations depends strongly on the autocorrela- 
tion time delay, which itself depends on the difference in 
arrival times at xi and X2 from scatterer s. For source lo- 
cations away from a straight line passing through xi and 
X2 (the "receiver line"), this difference in arrival time 
varies rapidly with source location causing fast varia- 
tions in the phase of the autocorrelations such that con- 
tributions from these source locations interfere destruc- 
tively. However, for source locations on the receiver line 
and not between xi and X2, the difference in arrival times 
is nearly constant regardless of the absolute location of 
s. The relative phase of these individual autocorrelations 
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is stationary, and the signals from these source locations 
interfere constructively. The S s (t) emanating from lo- 
cations on or very near the receiver line make the vast 
contribution to the cross correlation C(Xi,X2,?). 

From the physical, stationary phase argument above, 
one can already see that C(xi,X2,f) will contain infor- 
mation about G(xi,X2,*). Consider a source point s on 
the receiver line very near X2 but on the opposite side 
of X2 relative to Xi. Viewed from Xi, we could hardly 
distinguish whether a delta-function pulse from s origi- 
nated from s or from X2. The signal detected at Xi would 
be nearly equal to that caused by delta-function pulse 
from X2, which is exactly G(xi,x 2 ,f). A similar argu- 
ment can be made for the contributions from the other 
source terms on or near the receiver line but not neces- 
sarily close to X2. 

First applying a Fourier transform, Snieder lfl4ll works 
out the details of the stationary phase calculation by con- 
verting the sum over individual source terms s to an in- 
tegral over a source area density n(x). The result, given 
by his Eq. 24 and adapted to our setting, is 
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C(xi,x 2 ,(0) = — |S(co)| 2 [G(x 2 ,xi,co) / n(x)dx 



/CO 



G(xi,X2,co) / n(x)dx] 



(2) 



where C(xi,X2,G)) and G(xi,X2,Co) are the Fourier 
transforms of the cross correlation and Green's func- 
tion, |S(co)| 2 is power spectral density averaged over the 
sources S s (t), and c is the phase speed of EM-wave be- 
tween xi and X2. The integrals in Eq. [2] are carried out 
along two sections of the receiver line with Xi at x = 
and X2 at x — R. Snieder[14] performed his derivation 
for surface waves on a solid, and he accounted for differ- 
ent polarizations and different wave modes. EM waves 
are also two-dimensional waves, however, they are only 
scalar waves with a single mode. Therefore, in Eq.[2] we 
have dropped the polarization, mode indexes, and com- 
plex conjugation in Snieder's Eq. 24. Equation|2]is the 
starting point for our estimation of EM-wave Green's 
function from ambient frequency noise on an electrical 
interconnection. 



5 EM-wave Green's function esti- 
mation 

The electrical grid has been approximated as a two- 
dimensional continuum of transmission and generation, 
and frequency disturbances have been shown to prop- 
agate as EM waves as described above[3|. To use the 
Green's function estimation technique as described by 
Snieder lfHl . we must be careful in its application. In 



the continuum description of the grid, there are not em- 
bedded scatterers as described by Snieder| 14|, however, 
his scatterers simply acted as sources of uncorrelated in- 
jections of signals S s (t) from randomly placed locations. 
In the case of the electrical grid, random fluctuations of 
load at substation buses will play the same role. 

On short time scales, changes in load P at a bus in 
the transmission system are abrupt, and an individual 
change can be modeled as a step function in time. Previ- 
ously, we argued that such changes in load create instan- 
taneous frequency deviations at this bus (located at x s ) in 
the form of a Dirac delta-function in time, i.e. Af(x,t) = 
A/ V 8(x — x i )8(f). If we used the cross correlation in 
Eq.frito estimate G A /(xi,X2,f), the source terms S s (t) ~ 
Af s o[x — x s )8(f ) have a time-average component violat- 
ing a basic assumption made by Snieder[ 14]. To circum- 
vent this difficulty, we choose instead to work with the 
time derivative of the frequency deviation A/' and es- 
timate the Green's function G A ^/(xi,x 2 ,f). The source 
terms are S s (t) ~ Af s b(x — x s )db(t)/dt, which have no 
time-average component. The average of the power 
spectra in Eq. [5]is then |S(co)| 2 CO 2 . Substituting into 
Eq. [2] we find 



C A/ /(xi,x 2 ,co) 



/CO 



= -7tc[G A/ /(x 2 ,xi,co)A_ 

+ G A/ /(xi,x 2 ,«))A + ]. (3) 



Taking the inverse Fourier transform of Eq.[3]yields 

J C A f(x u x 2 ,t)dt = -7ic[G A //(x 2 ,xi,f)A_ 

+ G A/ /(x 1 ,x 2 ,f)A+], (4) 

where A_ and A + are the integrals over the source 
density in Eq. [2] from — °° — > and from R — > °o, re- 
spectively. If the electrical grid is in a quasi-steady 
state over the time window T, Af and Af satisfy 
the same homogenous partial differential equation and 
G A y/(xi,x 2 ,f) = G A /(xi,X2,f). Therefore, the time inte- 
gral of the cross correlation C A ^/(xi,x 2 ,f) in Eq. pflalso 
yields the Green's function for A/, a quantity of funda- 
mental importance. 

If we had a map of G A /(xi,x 2 ,f) for many pairs 
of locations throughout the interconnection, we could 
quickly compute the time dependence of the frequency 
deviation following a major disturbance to the grid's bal- 
ance of generation and load. In addition, G A /(xi,X2,f) 
can also be used to estimate the local EM-wave speed 
and attenuation throughout the interconnection which 
can then be used to quickly compute the frequencies and 
mode shapes of many inter-area oscillations and the full 
transient response to a grid disturbance. 
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Figure 1: The approximate location of the FNET sen- 
sors used in this analysis are indicated by the filled cir- 
cles. The black lines are the receiver lines for each pair 
of buses. The EM-wave speed for the receiver lines 
are computed from the arrival times estimated from the 
Green's function G\f(x\,X2,t) for each pair of buses. 

6 Source data and analysis 

To test the concepts described above, we have analyzed 
a small sample of data from the FNET WAMSQ. The 
FNET sensors provide time-synchronized frequency 
measurements every 0. 1 sees at many transmission buses 
throughout the U.S. Here, we focus on three closely- 
spaced buses in the Eastern interconnection shown in 
Fig. 1. The frequency measurements are actually per- 
formed at a location on the distribution system attached 
to the transmission bus, however, it has been shown 
that frequency disturbances travel through the distribu- 
tion system quickly compared to the transmission sys- 
tem making the FNET frequency measurements accu- 
rate representations of the bus-local frequency on the 
transmission bus|7|. 

Figure [2^ shows a 30-minute segment of frequency 
data from bus in Fig. 1 . On this time scale, the data 
from the other two buses closely track bus and are in- 
distinguishable. The data show normal variations about 
the base frequency, and important for our discussion, the 
data does not show any significant transients that can be 
attributed to a system disturbance. The data sets do show 
spurious spikes which we attribute to sensor noise. We 



have not made any specific effort to remove such noise. 

To apply Eq.|4] we must cross correlate dAf/dt at the 
three different sites. We compute the derivative at time f 
from the linear coefficient of a parabolic fit to a window 
of data Af wide and centered on f . The window At is var- 
ied in the following analysis to investigate the sensitivity 
of the results to this implicit filtering. Figure [2Jd shows 
the resulting dAf/dt for buses and 1 for Af ranging 
from 0.5 sees to 1.5 sees. We will use these two buses 
to investigate the how quickly the Green's function esti- 
mation technique converges. 

The cross correlation in Eq. [4] can be computed over 
a range of averaging times T. As T becomes longer, the 
cross correlation time integral averages over a greater 
number of load fluctuations and the destructive inter- 
ference between the S s (t) away from receiver lines be- 
comes more complete! 14], improving the estimation of 
GA/(xi,X2,f). However, a long T would limit our abil- 
ity to detect rapid changes in grid conditions in real 
time. Figure [3] shows the evolution in C^f{x\,X2,t) for 
T ranging from 1 minute to 16 minutes and for a At of 
1.0 sec and 1.5 sec for buses and 1. In Figure [5] T in- 
creases by a factor of two for each of the five curves in 
the groups of Af . 

Many FNET observations of the of the Eastern in- 
terconnection following a major disturbance show that 
the post-disturbance frequency is relatively constant for 
times greater than 10 sec after the disturbance!!*?! [lOl 171 . 
Since the cross correlation in Fig. [3] is the derivative of 
this response (i.e. of the Green's function), the cross 
correlation should approach a steady value of zero be- 
yond about 10 sec. For Af = 1.5 sec, the fluctuations in 
Coi beyond 10 sec are greatly diminished for T > 4 min 
and Coi appears to reach a relatively constant functional 
form. For Af = 1 .0 sec, T > 8 min is required to achieve 
similar behavior. These minimum values of T necessary 
to achieve adequate destructive interference of the noise 
sources off of the receiver line should be compared to 
the ~ 20 min averaging time to extract accurate results 
for bus-specific amplitudes of a single inter-area mode 
in JT3l . In principle, our method should require less av- 
eraging time because we are using all of the information 
in the cross correlation below a frequency of ~ 1/Af as 
opposed to the single frequency in fl3l . 

From the definition of the cross correlation in Eq. [T] 
and the interpretation in terms of Green's functions, the 
cross correlation in Fig.[3]should non-zero for both pos- 
itive and negative times shifts f, however, it should be 
zero at f = as a disturbance does not propagate with 
infinite speed. In Fig. [3] we only show the cross corre- 
lation for one half of the time axis, and it does not go 
to zero at f = 0. Both of these points can be explained 
by the distribution of sources around the receiver lines 
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Figure 2: a) A 30-minute window of frequency data 
from the FNET sensor near bus 0. The data shows nor- 
mal variability about the base frequency with no sig- 
nificant sudden increase or decrease in frequency. The 
downward spikes near 50 and 400 sec are likely just spu- 
rious, non-Gaussian sensor noise. The frequency data 
from buses 1 and 2 are indistinguishable from bus on 
this time scale, b) df/dt for bus (black) and 1 (red) for 
the first 20 sees of the window in a), df/dt vs t is calcu- 
lated from the linear coefficient a quadratic fit to the data 
in a) using different length time windows At centered on 
t. 



in Fig.[T]and the filtering we perform on df/dt. 

When a pulse is emitted at Xi, it propagates at finite 
speed and generates a response at X2 some time later, 
therefore, EM-wave propagation from xi — > X2 gener- 
ates a non-zero cross correlation in Eq. [TJ for t > 
and is due to sources contained in the integral A_ in 
Eq. |4] A similar argument concludes that EM-wave 
propagation from X2 —> X\ generates a non-zero cross 
correlation for t < and is due to sources contained in 
the integral A + . From these arguments, it is clear that 
Ga/'(x2,xi,?) is the causal Green's function for EM- 
wave propagation from Xi —> X2 and is non-zero for 
t > 0, and Ga/'(xi,X2,?) is the anti-causal Green's func- 
tion and is non-zero for t < 0. For time-reversal invari- 
ant systems, G&f(x2,X\,t) = Ga/'(xi,X2,— t), however, 
the cross correlation in Eq.|4]need not be symmetric be- 
cause the location and orientation of the receiver line in 
the network may yield weights A_ and A + that are not 
the same. Such a situation arises when one end of the 
receiver is near the edge of the network and proximity 
of the edge restricts either A_ or A + . 

For the three pairs of observation points in Fig.[T| one 
end of the receiver line is restricted in each case, i.e. 
near for — s- 1 and — !• 2 and near 2 for 1 — !• 2. For 
— >• 1 and — > 2, we use the cross correlation for t < 
and the anti-causal Green's function. For 1 — > 2, we use 
t > and the causal Green's function. These choices 
corresponds to waves emanating from noise source on 
the unrestricted end of the receiver line. It is these cross 
correlations (and resulting Green's functions) that are 
presented in Figs.[3]and[4] For waves from noise sources 
on the restricted end, we find cross correlations that are 
smaller than and lack the coherent oscillations of those 
in Figs. [3] and [4] 

The approximate symmetry of the cross correlation 
about t = also explains why the curves in Fig. [3] and 
|4]do not approach zero at t = 0. To obtain a reasonable 
signal-to-noise ratio, we filter df/dt using parabolic fits 
centered on t. The effect is to spread out any rapid time 
variation of df/dt causing sharp pulses to appear to ar- 
rive before causality would allow. Since our filtering 
window At is of the same order as the propagation time 
delay, a pulse is that is in reality sharp generates a con- 
tribution at t = when first passed through our filter. 

With a better understanding of the techniques for 
estimating the Green's function for EM wave propa- 
gation over an electrical grid, we compute the cross 
correlations for each pair of observation points using 
At = 1.5 sec and T = 4 min. The cross correlations are 
shown in Fig.|4^ and the integrals, i.e. the Green's func- 
tion GA/(xi,X2,f) for the local frequency, are shown in 
Figs. [4^ and b. The absolute scaling of C(xi,X2,f) and 
Ga/-(xi,X2,?) is arbitrary because we currently have no 
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Figure 3: The cross correlation of df jdt at buses and 
1, Cqi- The two groups of curves are computed using 
df/dt data generated using Af fitting windows 1.0 and 
1 .5 sec wide. Within each group, averaging times T of 
1, 2, 4, 8, and 16 minutes are used. The reduction in 
fluctuations of Coi beyond about 10 sec is used as an in- 
dication of the quality of destructive interference of the 
EM waves launched by sources away from the receiver 
line for buses and 1. These fluctuations become small 
and the functional form of Coi becomes relatively con- 
stant at Af = 1 .0 sec and T = 8 minutes or At — 1 .5 sec 
and T = 4 minutes. The zeroes for the Af = 1 .0 sec data 
are incrementally offset by 0. 1 along the vertical axis for 
clarity of the figure. The zeroes for Af = 1 .5 sec data are 
offset in a similar way with the zero for the T — I minute 
data at 1.0. 

method to estimate A + or A_ in Eq. [4] In an attempt 
to maintain the relative scaling between the three differ- 
ent observation pairs, we have used the product of root- 
mean-squares of df/dt at xi and X2 as estimates of A + 
and A_ 

As discussed earlier, GA/(xi,X2,f) is the frequency 
response at xi after a Dirac delta-function perturbation 
to the frequency at X2, which is equivalent to a sudden 
increase (decrease) in generation (load) at X2. Figure 
actually presents — Ga/(xi,X2,?) which is the response 
to a sudden loss of generation at X2. The general shape 
and time extent of the transient closely resembles ac- 
tual measurements (made with FNET) of such frequency 
declines after a loss of generation(9] [TO] [7J, however, 
the response in Fig. was estimated from ambient fre- 
quency noise. This noise can be analyzed in real time 
to provide an on-line prediction of the local and system- 
wide impact of such a major disturbance. 

The arrival time of an EM wave following a major 
disturbance has been estimated by finding the time of 
fastest frequency decline! 10 1 . For our analysis, this 
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Figure 4: a) The cross correlations for each pair of ob- 
servation points in Fig. [T] computed using Af = 1 .5 sec 
and T = 4 mintues. b) The negative of the integral of 
the cross correlations in a), which are proportional to 
— GA/(xi,X2,f), i.e. the frequency response at xi due to 
a sudden loss of generation at X2. 



is given by the time of the peak value of C(xi,X2,f). 
These arrival times are fo^i = 0.8 sec, fo^2 = 0.6 sec, 
and fi^2 = 0.4 sec. The straight-line geographic dis- 
tances between these buses are Lq^\ — 585 km, Lo^2 
160 km, and L\^2 = 510 km. The EM-wave speeds esti- 
mated from these values are shown in Fig.[TJ The slower 
EM-wave speeds appear to be correlated with regions of 
dense transmission and generation which is qualitatively 
consistent with the interpretation in [3 1. However, these 
analysis should be repeated for time windows T imme- 
diately following a major system disturbance, such as 
those in ifTOl . so that a quantitative comparison of the 
full GA/-(xi,X2,f) can be made. 
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7 Conclusions and future work 



References 



Using signal processing techniques developed for pas- 
sive imaging in seismology |fT4l . we have demonstrated 
how the transient response of the electrical grid to sud- 
den changes in load or generation can be estimated from 
ambient frequency noise gathered over a WAMS. Our 
method provides the entire time response, not just the 
phases and relative amplitudes for an individual EM 
standing-wave mode. We have shown how the quality 
of the Green's function for the transient response varies 
versus different levels of filtering of the raw signal and 
different cross-correlation integration times. In this ini- 
tial work, raw signal filtering with a time constant of 
about 1 — 1.5 sec is required to achieve an integration 
time of less than about 5 — 10 mintues. Shorter cross- 
correlation integration times are preferred as they allow 
for closer to real-time monitoring. 

There is much future work to do in this new area. 
Here, we mention a few possible directions: 

• PMU data has faster time resolution and different 
signal-to-noise characteristics. This method should 
be applied to PMU data to investigate if improve- 
ments in Green's function estimation are possible 

• The estimated Green's function should be com- 
pared to the transient responses generated by ma- 
jor disturbances such as loss of hundreds of MW 
of generation. The Green's function should be es- 
timated from the ambient frequency noise both be- 
fore and after the event and compared with the fre- 
quency decline during the event. 

• Using data from the FNET WAMS or PMU data, 
all nearest-neighbor Green's functions should be 
estimated and an EM-wave propagation speed and 
attenuation map developed for an entire intercon- 
nect. The mode shapes, frequencies, and damp- 
ing of all inter-area modes should be estimated and 
compared with off-line simulations. 
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